Phase reduction of stochastic Hmit cycle oscihators 
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We point out that the phase reduction of stochastic hmit cycle oscillators has been done incorrectly 
in the literature. We present a correct phase reduction method for oscillators driven by weak external 
white Gaussian noises. Numerical evidence demonstrates that the present phase equation properly 
approximates the dynamics of the original full oscillator system. 
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Many physical systems can be mathematically mod- 
eled by limit cycle oscillators. It is well known that the 
oscillator systems could exhibit a variety of behaviors. A 
fundamental theoretical technique for studying the oscil- 
lator dynamics is the phase reduction method (see e.g. 

This method has been widely and successfully ap- 
plied to coupled oscillators or an oscillator subjected to a 
regular external signal such as a periodic one. Consider- 
able theoretical progress has been made in understanding 
their dynamics by using this method. 

Recently, the dynamics of oscillators subjected to ex- 
ternal stochastic signals has also attracted much inter- 
est in connection with entrainment of independent os- 
cillators subjected to a common external noise. This 
common-noise-induced entrainment has been experimen- 
tally found in several systems as diverse as neuronal net- 
works [2] , ecological systems [1| , and lasers Q ■ Limit cy- 
cle oscillators driven by white Gaussian noise have been 
used as simple models for theoretically studying this en- 
trainment [1, S 0i 9- In these theoretical studies, the 
phase reduction method is applied to the noise-driven 
oscillators to derive a one dimensional equation for the 
phase variable only. Based on this phase equation, sev- 
eral reasonable theoretical results have been obtained. 
However, as we will show, the phase equation used in the 
above references is incorrect in the sense that in general 
it does not correctly describe the dynamics of the original 
full oscillator system even in the weak noise limit. 

The phase reduction is a powerful method for describ- 
ing the essential dynamics of oscillators. Application field 
of this method is expected to grow also in the case of 
stochastic oscillators. Therefore, it is essential to de- 
velop a phase reduction method for stochastic oscilla- 
tors. In the present paper, we consider a general class of 
limit cycle oscillators, which are subjected to white Gaus- 
sian noises, and develop the phase reduction method for 
these systems. Based on some numerical results, it is 
demonstrated that the present phase equation properly 
approximate the dynamics of the original full oscillator 
system while the incorrect version of phase equation fails. 
Finally, we make remarks on the results concerning the 
common- noise- induced entrainment, which have been ob- 
tained in Refs. [i,S0,[ii- 

Let X — {xi, . . . ,xn) & be a state variable vector 



and consider the equation 



X ^FiX) + G{X)at), 



(1) 



where F is an unperturbed vector field, G E R is a 
vector function, and ^(t) is the white Gaussian noise 
such that {^{t)) = and (C(i)C(s)) = 2DS{t - s), where 
(• • •) denotes averaging over the realizations of ^ and 6 
is Dirac's delta function. We call the constant D > 
the noise intensity. The noise-free unperturbed system 
X = F{X) is assumed to have a limit cycle with a fre- 
quency uj. We employ the Stratonovich interpretation for 
the stochastic differential equation ([1]). This interpreta- 
tion allows us to use the conventional variable transfor- 
mations in differential equations. A more general form 
of the noise term has been assumed in d, d, Q- However, 
we assume the form of Eq. ([1]) for simplicity. An exten- 
sion of the present derivation to a more general case is 
straightforward. 

Consider the unperturbed system X = and 
let Xq be its limit cycle solution. A phase coordi- 
nate (j) can be defined in a neighbourhood U of the 
limit cycle Xq in phase space. According to a conven- 
tional definition, we define the phase variable (f) so that 
(gradj(^0) • F{X) — uj may hold for any points in U. We 
can define the other — 1 coordinates r = (ri, . . . , r^-i) 
such that detd{(j),r)/dX ^ in U. We assume that 
r = a on the limit cycle, where a — (ai, . . . , a^r-i) 
is a constant vector. If we perform the transformation 
(xi, . . . , xn) (0, ri, . . . , rjv_i) in Eq. ([T]), we have the 
equation of the form 



^ = w + /i(</.,r)C(t), 

= M(l),r) + g,((j),r)^{t), 



(2) 
(3) 



where i = 1, . . . , N — 1. The functions h, fi, and gi are 
defined as follows: 

/i(0,r) = (grad,0).G(X(0,r)), (4) 
/,(</., r) = (grader,) •-F(X(0,r)), (5) 
5,(0, r) = (grader,) •G(X(0,r)), (6) 

where the gradients are evaluated at the point X{<j),r). 
These functions are periodic with respect to (fi: i.e., h{<f>-\- 
2ii,r) = h{(j),r), /i(0-h27r,r) = M(j),r), and g,(0 -t- 
27r,r) = g,{(j),r). 
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Equations ^ and ^ are Stratonovich stochastic dif- 
ferential equations. They can be converted into equiva- 
lent Ito stochastic differential equations. The cf) compo- 
nent of this Ito type equation is obtained as follows: 



D 



dh{(j), r) 



1^ —IZ. 



i=l 



dr. 



(7) 



In the case of weak noise < I?/cj <C 1, the deviation of 
r from a is expected to be small. Thus, we can use the 
approximation r = a in Eq. ([7]). Using this approxima- 
tion, we arrive at 

^ = uj + D[ Z'{cp)Z(cp) + W(cp) ] + Z{<P)m, (8) 

where 2^(0) and W{<j}) are given by 

Z(0) = /i(0,a), (9) 



(10) 



respectively. Since h and gi are periodic functions, Z(0) 
and W{4)) are also periodic: i.e., Z{4) + 2tt) — Z{(l)) and 
W{(j) + 2tt) — W{(j)). We may conclude that the reduced 
phase equation for the noise-driven oscillator ([T]) is given 
by Eq. ([8]). The oscillator dynamics is often studied 
by assuming a phase model instead of multidimensional 
differential equations. We emphasize that an equation of 
the form dS]) has to be assumed in studying the dynamics 
of oscillators with white Gaussian noises. 

In the previous studies the authors assumed 

the Ito type reduced phase equation of the form 



4) = uo + DZ{(t>)Z'{4)) + Z{(t))^{t). 



(11) 



Comparison of the present phase equation ^ and Eq. 
dTTjl clearly shows that the term DW{(j)) is dropped in 
the previously used equation pT|) . This term is 0(D) and 
is of the same order as DZ{ip)Z' (ip). Thus, in general, 
equation (fTT|) does not correctly describe the essential 
dynamics of the oscillator even in the lowest order ap- 
proximation as will be demonstrated. In the exceptional 
case ~ 0, it gives reasonable results. 

The reduced phase equation dH) is useful to calculate 
statistical quantities, which characterize the dynamics of 
oscillators subjected to white Gaussian noises. Using Eq. 
JS]) , we derive analytical expressions for two fundamental 
statistical quantities, which are the steady probability 
distribution of the phase variable and the mean frequency 
51. We will compare these quantities obtained by using a 
two dimensional oscillator model with those obtained by 
using its reduced phase model. 

Let P{(t), t) be the time-dependent probability distribu- 
tion function for the phase (p. The stochastic differential 
equation ([5]) is equivalent to the Fokker-Planck equation 



BP 8 

-Qi = -^[W + DZ{<P)Z'{<P) + DW{<P)}P] 



(12) 



We consider Eq. p2|) over the interval (p e [0, 27r] 
and assume the periodic boundary condition P(0, t) — 
P{2tt, t). The steady solution P{<p) is obtained by assum- 
ing dP/ dt ~ in Eq. . If we construct an asymptotic 
solution for P{(p) in the power of e = D/uj, then up to 
the first order we can obtain 

PW = ^ + ^ [Z{<P)Z'{^) - W{<P) +W]+ 0(e2), 

_ _ (12) 

where W is defined by = (27r)~i j^'^ W{(t))d(p. 

The mean frequency H of the oscillator is defined by 



fl = lim — 

T-oo T 



(p{t)dt. 



(14) 



This can be calculated by replacing the time average (|T4)) 
with the ensemble average: i.e., Vl ~ {(p). In the Ito equa- 
tion, unlike in Stratonovich formulation, the correlation 
between (p and ^ vanishes. If we take the ensemble aver- 
age of Eq. ((HI) , then we have 



{4>) = u: + D{Z{cp)Z' {cp) + W{cP)), 



(15) 



where we used the fact (Z((/>)^(i)) = {Z{(P)){i{t)) = 0. 
For an arbitrary function A{(p), the ensemble average can 
be calculated by using the steady probability distribution 
P{(p): i.e., {A) ^ J^"" A{(P)P{(P)d(p. If we use Eq. (fH]). 
we can obtain Q up to the second order in e as follows: 



1+eW + e^ 



{ZZ'Y-W^ + W -^0(e-^), (16) 



where {ZZ'Y = {lix)-^ J^'' {Z{(j>)Z' {(j))yd<p and = 

{2tt)~^ J^^ {W{(p)}^d(p. Since the white Gaussian noise 
has no characteristic frequency, intuitively, one might ex- 
pect that the noise does not cause any change in the oscil- 
lator frequency. However, this is not the case. Equation 
(I16p clearly shows that an external white Gaussian noise 
changes the mean frequency f2 in a general class of oscil- 
lators. It depends on the sign of W whether fl increases 
or decreases as the noise intensity increases. 

Equations (|13p and (fT6|) show that the term W{(p) in 
Eq. ^ significantly affects both the steady probability 
distribution P{(p) and the mean frequency in the first 
order of e. In particular, as shown by Eq. (fT6|) . the 
first order frequency shift is determined only from W{(p). 
Therefore, it is crucially important to include the term 
W{(p) into the reduced phase equation as in Eq. ^ . 
It is clear that the previously used phase equation (fTTj) 
cannot give proper approximations for P{(p) and f2. 

In order to validate the above phase reduction method, 
we carried out numerical calculations for an example of 
noise-driven oscillator. We compare P{(p) and fl between 
the theoretical and numerical results. We consider the 
Stuart-Landau (SL) oscillator 

i = x-coy~{x^+y^){x-c2y)+G^£.{t), (17) 
y = cox + y~{x^+y^){c2X + y)+Gy£.{t), (18) 
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where cq and C2 are constants, G = (Gx, Gy) is a vector 
function of (x, y), and ^ is the white Gaussian noise with 
the properties {£,{t)) = and (C(i)C(s)) = 2D 5{t - s). 
The SL oscillator has the limit cycle solution Xq = 
(coscji, sincji), where the natural frequency lo is given 
by a; = Co — C2. If we define the coordinates {(p, r) by the 
transformation 



+ C2 In r) , y 



+ C2lnr), (19) 



then (j) gives the phase variable and the limit cycle is 
represented by r = 1. 

The functions Z{(f>) and VF(0) can be obtained 
from Eq. p9)l . As examples, we consider the fol- 
lowing two types of G. Gi = (1,0) and G2 — 
(x, 0). For the first example, Z{(f)) and W{(l)) are 
given by Z{(l)) = — (sin (/) + C2 cos (/)) and W((/)) = 
{(1 + C2)/2} sin2(/). For the second example, they are 
given by Z{(j)) = — cos (sin (/) + C2 cos 0) and VF(0) = 
C2 cos^ (j) {~ cos 20 + C2 sin 2(f). Approximations for P{(t)) 
and n can be obtained by substituting these expressions 
for Z{(j)) and VF(0) into Eqs. ^ and As for 17, 

we can obtain {l/uj = 1 + 0{£^) for the first example Gi 
and n/uj = 1- (c2/4)£ + {(l + c^)/32}£2 + 0(£3) foj. ^he 

second example G2 ■ The former indicates that is inde- 
pendent of C2 and constant up to the second order in the 
first example. In contrast, the latter indicates that can 
either increase or decrease in the first order depending on 
C2 in the second example. 

In Figs. [Ha)-(d), numerical and theoretical results for 
P(0) are compared: filled circle and solid line represent 
P{4>) obtained by numerically solving Eqs. pT|) and pB]) 
and that given by Eq. p3|) . respectively. Theoretical 
predictions made by Eq. (jlip . which are obtained just by 
setting = in Eq. (fT5)) . are also shown by dashed line. 
Figures [1]^ a) and (b) are for the case of Gi while figures 
[Ijc) and (d) are for the case of G2. It is clear that the 
present phase model ([5]) gives precise approximations in 
all the cases. The agreements are excellent. In contrast, 
the incorrect version of phase equation pT|) does not give 
proper approximations at all in spite of the weak noise 
intensity. 

Figures [DJa) and (b) show the mean frequency plot- 
ted as a function of e = D/ui for the cases of Gi and 
G2, respectively. In all the numerical calculations, the 
natural frequency is set as m — \. The numerical results 
obtained by solving Eqs. P?|) and are shown by 
filled or open circle. The theoretical estimations given 
by Eq. (|16p are also shown by solid or dashed line. The 
theoretical estimation is Vl/uj — l + 0(e^) for Gi, which is 
independent of e and constant up to the second order in 
e. In Fig. UKa), the numerically obtained il is almost con- 
stant for (co,C2) = (1,0). This coincides with the above 
theoretical estimation. In the case of (co,C2) = (2,1), 
there is a deviation between the numerical and theoret- 
ical results: the numerical result shows an increase with 
increasing e. However, this increase in not linear with 
respect to e but a higher order one as shown in the in- 
set. In this sense, an agreement between the numerical 



(a) 






2345 



FIG. 1: Steady probability distribution P{(j>) of phase for 
noise-driven SL oscillator. Numerical result (•), analytical 
result Eq. (|13|l (solid line), and that obtained from Eq. pip 
(dashed line) are shown for e = 0.03. (a) Gi and (co,C2) = 
(1, 0); (b) Gi and (co, ca) = (2, 1); (c) G2 and (co, ca) = (2, 1); 
(d) G2 and (co,C2) = (0,-1). 




(b) 
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FIG. 2: Mean frequency Q. vs. e for (a) Gi and (b) G2. 
Numerical result (symbol) and analytical result Eq. (|16|) 
(line) are shown, (a) (co,C2) = (1,0) (•, dashed line) and 
(co,C2) = (2,1) (o, dashed line); (b) (co,C2) = (2,1) (•, 
dashed line) and (co,C2) = (0,-1) (o, solid line). The in- 
set in (a) is logarithmic plot of f2 — lj vs. e, where reference 
line for the scaling law is also shown. 



and theoretical results is confirmed up to the first or- 
der. In the case of G2, the theoretical estimation is given 
by n/uj ^ 1 - (c2/4)e + {(1 -f- c^)/32}e2 + 0{e^), which 
has a non- vanishing term of 0(e) except for C2 = 0. In 
Fig. mb), a good agreement between the numerical re- 
sult and this estimation is obtained in each of the cases 
(co, C2) = (2, 1) and (0, -1). If we use Eq. (fTT|) instead of 
Eq. ([8]), then we obtain the estimation fi/cj ~ l-|-0(£^), 
in which the 0{e) term vanishes. This estimation appar- 
ently disagrees with the numerical results. 

Figures. [T] and [2] clearly demonstrate that the reduced 
phase equation ([5]) precisely approximate the dynamics 
of stochastic oscillators with weak white Gaussian noises. 
In addition, it is apparent that the previously used equa- 
tion pT|) is erroneous. 

It is known that an entrainment could occur between 
two independent oscillators subjected to a common ex- 
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ternal white Gaussian noise. We discuss this entrainment 
phenomenon, paying particular attention to an effect of 
the term W{(j>). Consider the two equations 

X,^F{X,) + SF,{X,) + G,{Xi)^{t), 1^1,2, (20) 

where Xi G R^, F is an unperturbed vector field, SFi 
and SF2 are small deviations from it, Gi and G2 are 
slightly different vector functions, and £,{t) is the common 
white Gaussian noise with {^{t)) = and {^{t)^{s)) = 
2DS{t — s). The phase is defined by the unperturbed 
system X ~ F{X). Equation (PO]) can be reduced into 
the phase equation 

0, = + D[Zli^,)Z,i^,) + W^icj),)] + (21) 

where tOi represents the natural frequency. For simplicity, 
we assume uji is a constant. We introduce the average 
Z{(j)) = (Zi(0) + Z2(0))/2 and the difference (5Z(0) = 
(Zi(0) - Z2(0))/2. The difference SZ{(j)) is smah since 
Gi ~ G2 is assumed. 

Let 9 and ■0 be defined hy 6 = — (j>2 and ip = (j>i + 
^2 ~2u;t, where u is the average natural frequency defined 
by a; = {uji + u!2)/2. The variable 9 measures the phase 
difference between the two oscillators. In the case of weak 
noise, 9 and ip can be regarded as slow variables and 
thus the averaging approximation can be applied. If we 
perform the time-averaging and neglect the terms which 
are of the order of D|i5Z|, then we can obtain the Fokker- 
Planck equation for the probability distribution Q(t, 9, ip) 
as follows: 

^ = -{Su: + D{W,-W2)}^-D{W,+W2)^ 

+D-^He)Q]+D^[v{9)Q], (22) 

where Su = uji — uj2- The functions u and v are defined 
as u{9) = 2{r{0) - r{9)} and v{9) = 2{r(0) + r{9)}, 

where r{9) = i2n)-^ f^'' Z{<j))Z{<j)+9)d<j). Equation ^ 
has a steady solution Q{9), which is a function of 9 only. 
The entrainment phenomenon is characterized by Q{9). 
The steady solution is determined by the equation 

D^^^^^~{5u + D{W,~W2)]Q^C, (23) 



where C is an integration constant. In Eq. p6p . it has 
been shown that the mean frequency 51 shifts by DW 
in the lowest order. Equation indicates that this 
frequency shift effect appears as the effective detuning 
D(Wi-W2). The profile of Q{9) depends on the coeffi- 
cient of Q in Eq. ([23]). It has been shown that Q{9) has 
peaks at the zero points of u{9) and these peaks become 
narrower and higher, which corresponds to better syn- 
chronization quality, as the ratio {5u} -\- D{Wi — W2)}/D 
between the coefficient of Q and D becomes smaller 
Therefore, the profile of Q{9) depends on the functional 
form of Wi{(p) in Eq. (HI]). It may be concluded that 
the contribution due to W{<j)) in Eq. ([5]) is not negligible 
in the common- noise-induced entrainment. In addition, 
equation (j23p suggests that the synchronization quality 
could be improved if Sui and D{Wi — W2) cancel with 
each other. 

The previous works [1, H, 0, [1] have assumed the case 
Gi — G2. In this particular case, Wi = W2 holds and 
thus the effective detuning D{Wi — W2) vanishes in Eq. 
(1^5)) . Therefore, the same equation for Q{9) can be ob- 
tained even if Wi{(j)) in Eq. (PT]) is not taken into account. 
Because of this fact, fortunately, an analysis based on the 
incorrect phase equation ([TT]) also leads to correct results. 

In conclusion, we have developed the phase reduction 
method for a general class of limit cycle oscillators sub- 
jected to white Gaussian noises. Applying the present re- 
duced phase equation, we derived analytical expressions 
for the steady probability distribution P{cj)) of phase and 
the mean frequency Q. It has been found that an external 
white Gaussian noise gives rise to a frequency shift. We 
showed that these analytical estimations of P(0) and 
are in good agreement with numerical results to demon- 
strate that the present phase equation properly approxi- 
mates the dynamics of the original full oscillator system. 
In addition, we pointed out that an effect due to the fre- 
quency shift emerges also in the common-noise-induced 
entrainment. 

The authors would like to thank the members of NTT 
Communication Science Laboratories for their continual 
encouragements . 
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